/*******************************************************************************
Purpose: Replication code for main tables of Ouimet, Simintzi, & Ye; 
		 The Impact of the Opioid Crisis on Firm Value and Investment
*******************************************************************************/
set more off
set type double
clear all
cap log close
macro drop _all

*** NEED TO CHANGE THIS TO CORRECT FILE PATH ***
*cd "" /* Set directory here */

* Set paths:
global datapath "Data"
global outpath "Output"


************************************
//Table 1: Summary Statistics
************************************
//Panel A - County-level variables
*counties used in final sample
use "$datapath/Estab_variables.dta", clear
keep FIPS
duplicates drop FIPS, force

*merge raw county data between 2007 and 2015
merge 1:m FIPS using "$datapath/County_variables.dta", keep(match)
drop if opapreslg5==.

foreach var in Populationlg5 Incomelg5{
	gen `var'th=`var'/1000
}


*summary statistics use all county obserbvations between 2002 and 2010
*note: we use first-difference in regressions
matrix sumsts = J(8,4,0)
matrix rownames sumsts = "opapreslg5" "opapreser10lg5" "Populationlg5" "Incomelg5" "PopWhiteRatelg5" "Pop2064Ratelg5" "Pop65Ratelg5" "NeopMortalityRatelg5"
matrix colnames sumsts = "N" "Mean" "Median" "StDev"
local counter = 1
local Varlist "opapreslg5 opapreser10lg5 Populationlg5th Incomelg5th PopWhiteRatelg5 Pop2064Ratelg5 Pop65Ratelg5 NeopMortalityRatelg5"
foreach var in `Varlist'{
	winsor `var', gen(`var'w) p(0.01)
	qui sum `var'w, detail
	matrix sumsts[`counter',1] = r(N)
	matrix sumsts[`counter',2] = r(mean)
	matrix sumsts[`counter',3] = r(p50)
	matrix sumsts[`counter',4] = r(sd)
	local counter = `counter'+1
}
matrix list sumsts

sum opapreslg5, d


//Panel B: Establishment-level variable
*establishment used in final sample
use "$datapath/Estab_variables.dta", clear
keep siteid
duplicates drop siteid, force

*merge raw data between 2007 and 2015
merge 1:m siteid using "$datapath/IT", keep(match) keepusing(emple reven it_budget pcs Year)

*variables
foreach var in it_budget{
	gen `var'_revth=`var'/reven/1000
	gen `var'_empth=`var'/emple/1000
	gen `var'th=`var'/1000
	foreach v in `var'th `var'_revth `var'_empth{
		winsor `v', gen(`v'w) p(0.01)
	}
}

foreach var in pcs{
	gen `var'_rev=`var'/reven
	gen `var'_emp=`var'/emple
	foreach v in `var' `var'_rev `var'_emp{
		winsor `v', gen(`v'w) p(0.01)
	}
}

foreach var in emple reven{
	winsor `var', gen(`var'w) p(0.01)
}

*summary statistics use all establishment obserbvations between 2007 and 2015
matrix sumsts = J(4,4,0)
matrix rownames sumsts = "reven" "emple"  "it_budget" "pcs"
matrix colnames sumsts = "N" "Mean" "Median" "StDev"
local counter = 1
local Varlist "reven emple it_budgetth pcs"
foreach var in `Varlist'{
	qui sum `var', detail
	matrix sumsts[`counter',1] = r(N)
	matrix sumsts[`counter',2] = r(mean)
	matrix sumsts[`counter',3] = r(p50)
	matrix sumsts[`counter',4] = r(sd)
	local counter = `counter'+1
}
matrix list sumsts


************************************
// Table 2: Opioids and Individual Employment
************************************

// Robust #1 - Drop 10% high opioid doctors
use "$datapath/individual0.dta", clear
egen provprobrank = xtile(provprob), nq(10) by(year)
drop if provprobrank==10
save "$datapath/individual1.dta", replace

// Robust #2 - Drop counties with highest Gini index (10% sample)
use "$datapath/individual0.dta", clear
xtile ginirank= Gini, nq(10)
drop if ginirank==10
save "$datapath/individual2.dta", replace

// Robust #3 - ER 10
use "$datapath/individual0.dta", clear
keep if ertop10 == 1
save "$datapath/individual3.dta", replace


//Regression

use "$datapath/individual0.dta", clear

//Fixed effects
egen gfe = group(FIPS year dx age sex)
local fe "gfe plantyp"

//OLS
reghdfe emp opioid, keepsin a(`fe') cluster(FIPS)
outreg2 using "$outpath/Table2", keep(opioid) append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel

//First stage
reghdfe opioid provprob, a(`fe') cluster(FIPS)
outreg2 using "$outpath/Table2", keep(provprob) append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)', F, `e(F)') excel

//Second stage
ivreghdfe emp (opioid=provprob), a(`fe') cluster(FIPS)
outreg2 using "$outpath/Table2", keep(opioid) append bdec(3) tdec(3) nor2 addstat(Weakid F, `e(widstat)', R2, `e(r2)') excel



forvalue g=1/3{
use "$datapath/individual`g'.dta", clear

//Fixed effects
egen gfe = group(FIPS year dx age sex)
local fe "gfe plantyp"

//First stage
reghdfe opioid provprob, a(`fe') cluster(FIPS)
outreg2 using "$outpath/Table2", keep(provprob) append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)', F, `e(F)') excel

//Second stage
ivreghdfe emp (opioid=provprob), a(`fe') cluster(FIPS)
outreg2 using "$outpath/Table2", keep(opioid) append bdec(3) tdec(3) nor2 addstat(Weakid F, `e(widstat)', R2, `e(r2)') excel

}



************************************
//	Table 3: Opioids and Establishment Growth
************************************
use "$datapath/Estab_variables.dta", clear

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw opapreslg5chgw, a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table3", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
	reghdfe ln`depv'chgw opapreslg5chgw  `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table3", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}


************************************
//	Table 4: Opioids and IT Investment
************************************
use "$datapath/Estab_variables.dta", clear

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "it_budget it_budget_rev it_budget_emp pcs pcs_rev pcs_emp"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw opapreslg5chgw  `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table4", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}


************************************
//	Table 5: Robustness: Establishment Fixed Effects
************************************
use "$datapath/Estab_variables.dta", clear

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev it_budget_emp pcs pcs_rev pcs_emp"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw opapreslg5chgw `Controls5', a(naics4year ent_idgyear siteid) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table5", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}

************************************
//	Table 6: Heterogeneous Effects
************************************

//Panel A. The Labor Channel

use "$datapath/Estab_variables.dta", clear

*merge industry replaceability as of 2000
merge m:1 naics4 using "$datapath/replaceability_ind2000", nogen keep(master match)

*group by replaceability
foreach var in Replace2000{
	xtile `var'rank=`var', nq(2)
	gen `var'High=1 if `var'rank==2
	replace `var'High=0 if `var'High==. & `var'rank!=.
	drop `var'rank
}
gen op_Replace2000 = opapreslg5chgw * Replace2000High

*regression
set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "it_budget it_budget_rev it_budget_emp pcs pcs_rev pcs_emp"

foreach depvar in `DEPVAR'{
	reghdfe ln`depvar'chgw opapreslg5chgw op_Replace2000 `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table6_A", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}



//Panel B. Financial Constraints
use "$datapath/Estab_variables.dta", clear

xtile firmsizerank=firmsize, nq(2)
gen firmsizeHigh=1 if firmsizerank==2
replace firmsizeHigh=0 if firmsizeHigh==. & firmsizerank!=.
drop firmsizerank
gen op_firmsizeHigh = opapreslg5chgw * firmsizeHigh

egen FIPS_year=group(FIPS Year)

* without Firm FE
set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "it_budget it_budget_rev it_budget_emp pcs pcs_rev pcs_emp"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw opapreslg5chgw op_firmsizeHigh firmsizeHigh `Controls5', a(naics4year) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table6_B", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}




************************************
// Table 7: Robustness: Local Economic Conditions
************************************

///Panel A County-level income 
use "$datapath/Income.dta", clear
keep if Year == 2007 | Year == 2010
sort FIPS Year
by FIPS: gen Incomechg = Income - Income[_n-1]
keep if Year == 2010
keep FIPS Incomechg

* merge main dataset and drop lowest quantile ranked by housing price change
merge 1:m FIPS using "$datapath/Estab_variables.dta", nogen keep(match)
foreach v in Incomechg{
	xtile `v'rank=`v', nq(2)
	gen `v'low = 1 if `v'rank==1
	replace `v'low = 0 if `v'low==.
	gen op_`v'low = opapreslg5chgw * `v'low
	local interaction "op_`v'low opapreslg5chgw `v'low"

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev pcs pcs_rev"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw `interaction' `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table7_A", addtext(Firm-Year FE, Yes, Ind-Year FE, Yes, low, below med, income, ch 2007-2010) append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}
}



///Panel B Housing price index
use "$datapath/HPI.dta", clear
keep if Year == 2007 | Year == 2010
sort FIPS Year
by FIPS: gen hpichg = hpi - hpi[_n-1]
keep if Year == 2010
keep FIPS hpichg

* merge main dataset and drop lowest quantile ranked by housing price change
merge 1:m FIPS using "$datapath/Estab_variables.dta", nogen keep(match)
foreach v in hpichg{
	xtile `v'rank=`v', nq(2)
	gen `v'low = 1 if `v'rank==1
	replace `v'low = 0 if `v'low==.
	gen op_`v'low = opapreslg5chgw * `v'low
	local interaction "op_`v'low opapreslg5chgw `v'low"


set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev pcs pcs_rev"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw `interaction' `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table7_B", addtext(Firm-Year FE, Yes, Ind-Year FE, Yes, low, below med, HPI, ch 2007-2010) append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}
}



///Panel C County employment
use "$datapath/Employment.dta", clear
keep if Year == 2007 | Year == 2010
sort FIPS Year
by FIPS: gen employedchg = employed - employed[_n-1]
keep if Year == 2010
keep FIPS employedchg

* merge main dataset and drop lowest quantile ranked by housing price change
merge 1:m FIPS using "$datapath/Estab_variables.dta", nogen keep(match)
foreach v in employedchg{
	xtile `v'rank=`v', nq(2)
	gen `v'low = 1 if `v'rank==1
	replace `v'low = 0 if `v'low==.
	gen op_`v'low = opapreslg5chgw * `v'low
	local interaction "op_`v'low opapreslg5chgw `v'low"

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev pcs pcs_rev"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw `interaction' `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table7_C", addtext(Firm-Year FE, Yes, Ind-Year FE, Yes, low, below med, emp, ch 2007-2010) append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}
}




************************************
//	Table 8: Robustness: Excluding Pill Mill Counties and Florida
************************************
//Panel A - exclude pill mill counties
use "$datapath/Estab_variables.dta", clear

merge m:1 FIPS using "$datapath/arcos_co_2006.dta", nogen keep(master match)

foreach var in co_mme_top5_count2006{
	xtile `var'rank=`var', nq(4)
}
drop if co_mme_top5_count2006rank==4

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev pcs pcs_rev"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw opapreslg5chgw `Controls5', a(naics4year ent_idgyear) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table8_A", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}



//Panel B - exclude Florida
use "$datapath/Estab_variables.dta", clear

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev pcs pcs_rev"

foreach depv in `DEPVAR'{
	reghdfe ln`depv'chgw opapreslg5chgw `Controls5' if STFIPS!="12", a(naics4year ent_idg) cluster(FIPS ent_idg) tol(0.0001)
	outreg2 using "$outpath/Table8_B", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}



************************************
//	Table 9: 2SLS: Emergency Room Opioids Instrument
************************************
use "$datapath/Estab_variables.dta", clear

set more off
local Controls5 "lnIncomelg5chgw lnPopulationlg5chgw PopWhiteRatelg5chgw Pop2064Ratelg5chgw Pop65Ratelg5chgw NeopMortalityRatelg5chgw"
local DEPVAR "reven emple it_budget it_budget_rev it_budget_emp pcs pcs_rev pcs_emp"

reghdfe opapreslg5chgw opapreser10lg5chgw `Controls5', a(ent_idgyear naics4year) cluster(FIPS) tol(0.0001)
outreg2 using "$outpath/Table9", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel

foreach depv in `DEPVAR'{
	ivreghdfe ln`depv'chgw `Controls5' (opapreslg5chgw = opapreser10lg5chgw), a(ent_idgyear naics4year) cluster(FIPS)
	outreg2 using "$outpath/Table9", append bdec(3) tdec(3) nor2 addstat(R2, `e(r2)') excel
}


************************************
// Table 10: Abnormal Returns around the Passage of State Laws on Opioids
************************************
* Note: Return data is not real.

use "$datapath/return_data.dta", clear

*regression
foreach v in retp retp_high_HQ_empl retp_low_HQ_empl retp_high_pc retp_low_pc{
	reghdfe `v'rf EventWindow mktrf smb hml, a(eventdate) cluster(eventdate)
	outreg2 using "$outpath/Table10", append bdec(5) tdec(5) nor2 addstat(Adjusted R2, `e(r2_a)') excel	
	reghdfe `v'rf EventWindow mktrf smb hml umd, a(eventdate) cluster(eventdate)
	outreg2 using "$outpath/Table10", append bdec(5) tdec(5) nor2 addstat(Adjusted R2, `e(r2_a)') excel	
}



************************************
//	Figure 1: Map of Opioid Prescriptions
* We use software Tableau to produce the map
************************************
use "$datapath/opioids.dta", clear
keep if Year >= 2001 & Year <= 2010
bysort FIPS: egen op_avg = mean(opapres)
keep FIPS op_avg
duplicates drop FIPS, force












